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ABSTRACT 

An  investigation  of  inertial  coupling  and  its  contribution  to  wing  rock  in  the 
F-14A  aircraft  has  been  conducted.  Wind  tunnel  data  was  used  to  obtain  the 
stability  parameters  for  angles  of  attack  from  zero  to  25  degrees,  after  which 
linear  and  nonlinear  analyses  of  the  equations  of  motion  were  completed.  The 
linearized  analysis  of  the  uncoupled  longitudinal  and  lateral-directional  equations 
was  included  to  provide  a  baseline  for  comparison  with  the  fully  coupled, 
nonlinear  equations.  In  both  cases,  the  equations  of  motion  were  solved 
numerically  and  time  history  traces  produced  to  illustrate  aircraft  response. 
Results  indicate  that  a  stable  short  period  mode  can  feed  damping  energy  into  an 
unstable  dutch  roll  mode  via  the  coupling  of  the  equations  to  produce  a  stable 
limit  cycle  very  similar  to  those  experienced  in  the  aircraft,  Numerous 
suggestions  for  follow  on  research  are  presented. 
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I.  INTRODUCTION 


Much  has  been  written  on  the  subject  of  linearized  aircraft  flight  mechanics. 
Under  the  assumption  of  small  perturbations,  the  equations  of  motion 
representing  rigid  body  aircraft  dynamics  can  be  reduced  to  two  sets  of 
uncoupled,  linear  differential  equations  with  constant  coefficients  which  can  be 
readily  solved  to  yield  characteristic  frequencies,  damping  ratios,  and 
amplitude/phase  relationships  of  the  aircraft's  "natural  modes",  as  well  as  time 
history  traces  of  aircraft  response  to  these  modes  due  either  to  initial  conditions 
or  control  inputs.  The  approximations  made  in  thi.  analysis  of  the  linearized 
equations  are  quite  good  as  long  as  the  aircraft  motion  is  analyzed  within  the 
boundaries  of  linear  behavior  (i.e.,  small  angles).  Under  the  assumptions  made 
in  this  relatively  benign  region,  the  longitudinal  aircraft  motions  are  isolated 
from  the  lateral-directional  motions.  Analytical  results  from  linearized  theory 
agree  well  with  flight  tests  in  this  regime,  with  negligible  coupling  occurring  in 
flight.  As  good  as  these  approximations  are,  they  break  down  completely  at  high 
angle  of  attack,  sideslip  angle  or  high  angular  rates  and  fail  to  provide  accurate 
results  when  analyzed  in  this  manner.  Thus,  the  full  set  of  nonlinear,  coupled 
equations  of  motion  must  be  analyzed  to  obtain  adequate  results.  Whereas  the 
motion  and  inertially  related  nonlinear  terms  were  neglected  in  small 
perturbation  theory,  these  terms  provide  coupling  between  the  longitudinal  and 
lateral-directional  dynamics  in  the  nonlinear  analysis,  altering  the  overall 
dynamic  response. 

The  differences  between  linear  and  nonfnear  flight  dynamics  described 
above  form  the  foundation  for  this  study.  In  particular,  this  analysis  probes  into 
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the  nonlinear  phenomena  known  as  wing  rock.  Numerous  studies  conducted 
over  recent  years  in  an  attempt  to  identify  the  specific  cause(s)  of  wing  rock  have 
highlighted  aerodynamic  hysteresis  [Ref.  1],  complex  asymmetric  vortex 
interactions  for  slender  deltas  [Ref.  2,  3]  and  inertial  cross  coupling  [Ref.  4]  as 
potential  candidates  for  the  wing  rock  motion.  The  validity  of  these  studies  is 
certainly  not  under  question;  however,  not  one  of  them  can  be  chosen  as  the 
specific  cause  for  wing  rock  in  any  aircraft  without  careful  consideration  due  to 
the  highly  configuration-dependent  nature  of  the  problem.  It  is  with  this  in  mind 
that  this  study  investigated  the  effects  of  inertial  cross  coupling  as  a  contributor 
to  the  wing  rock  motion  in  the  F-14A  aircraft.  Although  wing  rock  has  been 
reported  in  many  other  aircraft  (e.g.,  A-4,  F-4,  F-5,  [Ref.  1];  HP-115  [Ref.  2]; 
T-38  [Ref.  4];  F-15  [Ref.  5];  F/A-18,  X-29  [Ref.  6])  this  aircraft  was  modeled  for 
a  variety  of  reasons.  Perhaps  most  important  among  these  reasons  was  the 
rapidly  growing  emphasis  on  high  angle  of  attack  research  and  technology,  along 
with  its  potential  impact  on  the  future  combat  capability  of  the  F-14.  The  very 
existence  of  aircraft  currently  displaying  excellent  high  angle  of  attack 
performance  such  as  the  F-16,  F/A-18,  and  Mig-29  merely  suggest  that  future 
designs  will  display  even  greater  capability  in  this  regime.  As  the  F-14  enters  its 
third  decade  of  service,  it  carries  with  it  a  battle-proven  record  of  superior 
performance.  But  the  realities  of  reduced  military  spending  and  a  huge  federal 
deficit  make  it  all  too  clear  that  every  weapon  system  in  the  inventory  will  be 
used  to  the  very  end  of  its  service  life  (and  as  we  are  witnessing  with  A-6's,  at 
times  beyond  the  projected  service  life  by  sending  aircraft  through  rework 
facilities  for  structural  repair  and/or  enhancement;  a  less  time  consuming  and 
less  expensive  alternative  to  development,  test  and  acceptance  of  new  designs). 
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Thus,  if  trends  in  the  defense  contracting  industry  and  weapons  procurement 
branches  of  the  armed  services  continue  to  proceed  as  they  have  in  the  80's,  it  is 
not  unreasonable  to  assume  that  the  F-14  will  remain  in  service  for  the 
foreseeable  future,  meaning  at  least  into  the  2ist  century.  The  F-14  will 
undoubtedly  face  aircraft  with  superior  high  angle  of  attack  capability  and  may 
be  forced,  depending  on  rules  of  engagement,  current  tactics,  etc.,  to  fight  these 
aircraft  in  the  high  AOA  regime.  Even  with  the  current  upgrades  to  engines  and 
avionics,  the  basic  airframe  of  the  F-14  will  remain  the  same.  Therefore,  the 
dynamic  behavior  of  the  aircraft  from  a  stability  standpoint  will  also  remain 
unchanged  unless  alterations  are  made  to  the  aircraft's  stability  augmentation 
system.  If  an  F-14  experiences  wing  rock  while  engaged  in  an  aerial  encounter, 
especially  at  very  high  angle  of  attack,  the  only  current  recourse  is  to  neutralize 
lateral  and  directional  controls  and  momentarily  reduce  the  angle  of  attack. 
Obviously,  this  course  of  action  may  compromise  the  aircraft  and  prove  to  be 
unrecoverable  for  the  crew.  If  it  is  within  our  means  to  understand  this 
phenomenon,  then  we  are  that  much  closer  to  providing  and  implementing  a 
solution  in  the  form  of  control  laws  to  avoid  it. 

Another  important  reason  for  studying  the  F-14  was  the  relatively  simple 
flight  control  system  incorporated  into  its  design.  As  opposed  to  many  of  the 
newest  high  performance  aircraft,  the  F-14  flight  control  system  uses  only  pilot 
input  via  mechanical  linkages  and  a  three  axis  stability  augmentation  system 
(SAS)  to  control  the  aircraft  during  normal  maneuvering.  The  aircraft  can  be 
flown  freely  within  its  envelope  without  pitch  SAS  and  is  restricted  to  roll  SAS 
off  above  15  units  AOA  for  subsonic  maneuvering.  Thus,  only  the  function  of 
the  yaw  SAS,  the  most  critical  of  the  three,  was  of  concern  to  this  study. 
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Although  one  design  feature  of  the  yaw  SAS  was  to  provide  increased  directional 
stability  at  high  angle  of  attack,  the  aircraft  can  still  be  maneuvered  normally 
without  it  "if  extra  care  is  taken  to  control  yaw  excursions  with  rudder."  [Ref. 
7:p.  IV-11-1]  Therefore,  conducting  an  analysis  of  the  stick  fixed  stability 
characteristics  of  this  aircraft  without  knowledge  of  the  details  of  its  stability 
augmentation  system  is  not  unreasonable.  For  aircraft  with  "fly-by-wire"  flight 
control  systems,  the  results  of  a  similar  stability  analysis  may  vary  by  a  wider 
margin  (from  actual  aircraft  response)  due  to  the  influence  of  the  flight  control 
computer  on  the  control  surfaces  in  response  to  changing  flight  conditions,  even 
with  no  pilot  input. 

The  F-14  also  proved  to  be  an  excellent  choice  due  to  the  availability  of  an 
extensive  aerodynamic  data  base.  Lastly,  an  abundance  of  pilot  reports  from 
fellow  aviators  at  NFS  was  readily  available. 

In  an  effort  to  ultimately  obtain  numerically  derived  time  history  traces 
illustrating  wing  rock  for  the  F-14,  the  sections  that  follow  provide  a  brief 
background  on  F-14  high  AO  A  flight  characteristics,  the  aerodynamic  data  base 
from  which  the  aircraft  stability  derivatives  were  obtained  and  a  summary  of  the 
equations  of  motion.  The  study  continues  by  describing  the  numerical 
procedures  used,  then  highlights  the  linearized  and  nonlinearized  results.  A 
number  of  recommendations  for  related  research  are  presented  following  the 
conclusions. 
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n.  F-14A  FLIGHT  CHARACTERISTICS  REVIEW 

A.  NATOPS  FLIGHT  MANUAL  REVIEW 

The  F-14A  NATOPS  Flight  Manual  [Ref.  7]  mentions  wing  rock  and 
coupling  tendencies  in  the  Flight  Characteristics  section  of  the  manual.  The 
following  quotations  from  the  manual  are  provided  to  familiarize  the  reader  with 
the  nonlinear  dynamics  of  this  aircraft. 


Coupling  occurs  when  motions  in  more  than  one  axis  interact.  The  F- 
14A,  like  all  high-performance  airplanes  capable  of  producing  high  rated, 
multiple  axes  motion,  is  susceptible  to  coupling.  High  rate,  multiple  axes 
motions  particularly  at  high  AOA  can  produce  violent  coupled  departures. 
[Ref.7:p.  lY-ll-lS) 

Although  the  maneuver  slats  increase  the  severity  of  the  wing  rock 
between  20  and  28  units  AOA,  overall  departure  resistance  of  the  aircraft  is 
greatly  improved.  The  wing  rock  may  be  damped  with  rudders,  but  greater 
difficulty  may  be  encountered  with  maneuver  slats,  especially  at  low 
airspeeds.  If  this  occurs,  the  wing  rock  may  be  damped  by  neutralizing  the 
lateral  and  the  directional  controls  and  momentarily  reducing  AOA  to 
below  20  units.  [Ref.  7:p.  IV-11-71 

At  20  to  28  units  AOA,  reduced  directional  stability  is  apparent,  and 
even  small  control  inputs  will  cause  yaw  oscillations  that,  if  unchecked,  can 
produce  a  mild  wing  rock  (+/*  10  to  15  degrc^i  s)....  Wing  rock  at  20  to  28 
units  AOA  will  be  more  severe  < +/-  25  degrees)  and  more  difficult  to  damp 
with  maneuver  slats  extended  (due  to  incrersed  dihedral  effect).  [Ref.  7:p. 
IV-11-11] 
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In  the  takeoff  and  landing  configuration,  wing  rock  is  also  experienced  at 
high  AOA  during  approach  to  stall,  as  indicated  by  the  following  excerpt  from 
the  NATOPS. 

At  25  units  AOA  divergent  wing  rock  and  yaw  excursions  define  the 
stall.  Sideslip  angle  may  reach  25  degrees,  and  bank  angle  90  degrees 
within  6  seconds  if  AOA  is  not  lowered....  Extending  the  speed  brakes  ... 
improves  directional  stability  significantly,  reducing  the  wing  rock  and  yaw 
tendency  at  25  units  AOA.  Stall  approaches  should  not  be  continued  beyond 
.  the  first  indication  of  wing  rock.  When  wing  rock  occurs,  the  nose  should 
be  lowered  and  no  attempt  should  be  made  to  counter  the  wing  rock  with 
lateral  stick  or  rudder.  (Ref.  7:p.  IV-11-18] 

B.  FLIGHT  TEST  TIME  HISTORIES 

To  further  illustrate  the  characteristics  of  wing  rock  for  this  aircraft,  several 
time  history  traces  from  two  actual  flight  tests  are  shown  in  Figures  la- Id  and 
2a-2d  [Ref.  8].  In  both  sets,  it  can  be  noted  that  the  lateral  acceleration  of  the 
c.g.  shown  in  Figures  Id  and  2d  can  be  converted  into  a  sideslip  perturbation. 

The  first  set  of  traces.  Figures  la- Id,  was  obtained  from  a  flight  test 
conducted  at  approximately  20,000  ft.,  Mach  .32,  with  maneuver  flaps  extended, 
wings  swept  fully  forward  and  the  SAS  off.  The  initial  angle  of  attack  (ARI 
ALPHA)  from  the  plot  was  14  degrees.  This  angle  corresponded  to 
approximately  12  degrees  true  AOA  via  the  conversion 

AOAtrue  =  -8122*AOAARI  +  .7971  (01) 
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[Ref.  8]  and  17  units  indicated  AO  A  via  the  conversion 


AOAunits  =  1 .0989  *(AOAtrue  +  3.01 )  (02) 

which  is  valid  for  the  F-14  operating  at  Mach  number  less  than  .4  [Ref.  9]. 

As  the  time  history  began,  the  aircraft  was  decelerating  and  most  likely 
experiencing  light  buffet  (note  the  rapid  variations  in  normal  acceleration, 
Figure  la).  In  order  to  maintain  altitude,  the  pilot  commanded  aft  stick.  This 
input  resulted  in  corresponding  increases  in  stabilator  position,  angle  of  attack 
and  pitch  attitude  (Figure  lb).  Around  15  seconds  and  at  the  equivalent  of 
approximately  18.5  units  AOA,  very  slight  roll  rate  oscillations  began  which 
resulted  in  a  20  degree  left  wing  down  attitude  after  about  28  seconds  and 
coupled  into  small  angle  of  attack  and  sideslip  perturbations.  The  pilot 
responded  with  a  slight  right  stick  correction  towards  wings  level,  at  which  time 
larger  roll  rate  and  roll  amplitude  oscillations  began  to  develop,  at  around  30-35 
seconds  (Figure  Ic).  The  angle  of  attack  at  the  point  where  pronounced  wing 
rock  oscillations  first  appeared  was  the  equivalent  of  20.5  units.  At  the  same 
time,  larger  sideslip  and  yaw  rate  oscillations  began  to  appear  (Figure  Id),  as 
well  as  a  rapid  nose  down  pitch  rate.  As  the  left  wing  down  oscillations  were 
arrested  and  the  mean  value  of  the  wing  oscillation  returned  towards  wings  level, 
the  lateral  stick  input  was  taken  out.  With  the  decrease  in  pitch  attitude  came  a 
loss  in  altitude  and  an  increase  in  velocity.  Recovery  began  by  applying  forward 
stick  to  reduce  angle  of  attack  at  about  45  seconds,  although  the  oscillations 
continued  beyond  the  end  of  the  trace  at  60  seconds. 
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Figure  1  a.  Flight  Test  Time  History  Traces 
F-14A  IX.  BUNO  157991 ,  FLT  264,  GW  55386  LB.,  C.G.  6.8  %  MAC 
SWEEP=19  DEG.  FLAPS=1 1  DEG,  ALT=20200  FT.  MACH  .32,  SAS  OFF 
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Figure  lb.  Flight  Test  Time  History  Traces 
F-14A  IX.  BUNO  157991,  FLT  264,  GW  55386  LB.,  C.G.  6.8  %  MAC 
SWEEP=19  DEO,  FLAPS=1 1  DEO,  ALT=20200  FT,  MACH  .32,  SAS  OFF 
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Figure  Ic,  Flight  Test  Time  History  Traces 
F-14A  IX,  BUNO  157991,  FLT  264,  GW  55386  LB.,  C.G.  6.8  %  MAC 
SWEEP=19  DEG,  FLAPS=1 1  DEG,  ALT=20200  FT,  MACH  .32,  SAS  OFF 
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Figure  Id.  Flight  Test  Time  History  Traces 
F-14A  IX.  BUNO  157991,  FlT  264,  GW  55386  LB.,  C.G.  6.8  %  MAC 
SWEEP=19  DEG.  FLAPS=1 1  DEG,  ALT=20200  FT,  MACH  .32.  SAS  OFF 
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The  second  set  of  traces,  Figures  2a-2d,  was  taken  from  the  same  test  flight 
at  a  nearly  identical  flight  condition,  except  that  the  altitude  was  approximately 
19,000  ft.  The  initial  angle  of  attack  (ARI  ALPHA)  from  the  plot  was  17 
degrees.  This  corresponded  to  approximately  14.5  degrees  true  AOA  and  19.5 
units  as  seen  in  the  cockpit.  In  this  demonstration,  the  aircraft  began  by  quickly 
decelerating  and  increasing  AOA  to  19  degrees  ARI  (approximately  16.5  degrees 
true  and  21  units,  Figures  2a  and  2b).  The  immediate  presence  of  roll  and  roll 
rate  oscillations  without  lateral  control  inputs  is  evident  in  Figure  2c,  as  well  as 
sideslip  and  yaw  rate  oscillations  in  Figure  2d.  The  roll  amplitude  eventually 
reached  approximately  +/-  45  degrees  with  a  period  of  about  5  seconds.  It  is 
evident  that  frequency  doubling  in  the  AOA  penurbation  occurred,  as  the  period 
was  about  2*2.5  seconds.  The  AOA  reached  approximately  25.5  units  at  its 
highest  point.  Also,  a  pronounced  decrease  in  pitch  angle  was  evident  once  small 
AOA  perturbations  appeared  at  around  10  seconds.  Furthermore,  loss  of  altitude 
appeared  as  the  roll  angle  oscillations  built  beyond  +/-  30  degrees  from  the  mean 
value  at  around  25  seconds. 

All  of  these  illustrations  clearly  indicate  that  the  motions  which  occur  at  high 
angle  of  attack  are  quite  complex  and  that  nonlinear  coupling  occurs  between  the 
longitudinal  and  lateral-directional  equations  of  motion. 
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Figure  2a.  Flight  Test  Time  History  Traces 
F-14A  IX,  BUNO  157991,  FLT  264,  GW  55345  LB.,  C.G.  7.1  %  MAC 
SWEEP=19  DEG.  FLAPS=1 1  DEG.  ALT=18988  FT,  lACH  .32,  SAb  OFF 
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Figure  2b.  Flight  Test  Time  History  Traces 
F-14A  IX.  BUNO  157991,  FLT  264,  GW  55345  LB.,  C.G.  7.1  %  MAC 
SWEEP=19  DEG.  FLAPS=1 1  DEG,  ALT=18988  FT,  MACH  .32.  SAS  OFF 
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Figure  2c.  Flight  Test  Time  History  Traces 
F-14A  IX,  BUNO  157991,  FLT  264,  OW  55345  LB.,  C.O.  7.1  %  MAC 
SWEEP=19  DEO,  FLAPS=1 1  DEO,  ALT=1898f;  FI',  MACH  ,32,  SAS  OFF 
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Figure  2d,  Flight  Test  Time  History  Traces 
F-14A  IX,  BUNO  157991,  FLT  264,  OW  55345  LB..  CO.  7.1  %  MAC 
SWEEP=19  DEO,  FLAPS=11  DEO,  ALT=18988  FT,  MACH  .32,  SAS  OFF 


m.  A  DESCRIPTION  OF  THE  AERODYNAMIC  DATA  BASE 


The  stability  derivatives  needed  to  obtain  the  coefficients  in  the  equations  of 
motion  for  the  F-14  came  from  an  aerodynamic  data  base  provided  by  the 
Grumman  Corporation,  Bethpage,  New  York  [Ref.  9].  The  data  base  consisted 
of  aerodynamic  and  stability  coefficients  from  wind  tunnel  and  spin  tunnel  tests 
conducted  on  a  flight  test  aircraft  at  the  NASA  Langley  Research  Facility.  The 
tabulated  data  represented  two  different  flight  regimes;  a  low  speed  regime 
where  Mach  number  remained  at  or  below  M-.6  and  a  high  speed  regime  for 
Mach  number  greater  than  M=.6.  Typically,  the  low  speed  data  appeared  in  the 
tables  as  a  function  of  angle  of  attack  and  sideslip,  although  some  stability 
derivative  coefficients  were  dependent  solely  on  AO  A,  while  others  were 
dependent  upon  AOA,  sideslip  and  control  surface  deflection.  All  data  in  the 
tables  was  referenced  to  the  body  axis  coordinate  system  (x  •  axis  corresponding 
to  the  aircraft's  longitudinal  axis).  Therefore,  coordinate  transformations 
through  the  angle  of  attack  were  required  to  obtain  coefficients  related  to  lift  and 
drag.  Additionally,  crossplots,  curve  fits  and  a  number  of  other  techniques  were 
employed  to  obtain  all  of  the  necessary  information  for  the  stability  analysis.  A 
detailed  description  of  the  techniques  used  to  obtain  the  stability  parameters  from 
the  wind  tunnel  data  appears  in  Appendix  A. 

It  is  very  important  to  remember  the  distinction  between  true  angle  of  attack 
expressed  in  degrees  and  in  units.  The  conversion  between  degrees  true  and 
units  for  the  F-14  at  Mach  number  less  than  .4  was  given  in  equation  (02).  The 
net  effect  of  this  conversion  is  to  subtract  five  from  the  indicated  AOA  expressed 
in  units  to  get  AOA  in  degrees  (’.e.,  25  units  =  20  degrees).  This  distinction 
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should  be  kept  in  mind  when  comparing  the  NATOPS  flight  characteHstic5  to  the 
results  of  this  study,  which  indicate  aircraft  response  at  various  tnac  angles  of 
attack  expressed  in  degrees. 
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IV.  EQUATIONS  OF  MOTION  DEVELOPMENT 


The  equations  of  motion  describing  the  dynamic  behavior  of  an  aircraft  are 
typically  developed  in  a  rotating  reference  frame  Hxed  to  the  aircraft's  body 
4txes.  The  aircraft  is  assumed  to  be  a  rigid  body,  such  that  aeroelastic  and 
gyroscopic  effects  are  neglected.  Additionally,  the  variation  in  aircraft  velocity 
is  assumed  to  be  negligible  for  the  purposes  of  this  study.  This  approximation 
allows  for  the  removal  of  the  X  force  equation  and  the  perturbation  velocity 
from  the  longitudinal  equations  of  motion.  After  the  assumption  of  small 
perturbations  is  made  and  the  linearization  of  the  equations  is  completed,  they 
can  be  conveniently  expressed  in  state  space  format  as  shown  below. 

Xiong  =  Along  Xiong  (03) 

and 

Xlat-dir  =  A|at-dir  Xiat-dir  (04) 

where 

Xiong  =  [  Cl  q  0  J  (05) 
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Mg'  s  Mg  +  y  (Mq  Zg) 


Mq*  =  Mq  +  Mq 


Xiat-dir 


=  [  p  p  (^  r] 


(07) 

(08) 

(09) 


^lat-dir  = 


Yp  gcos80 

U  ®  U 


Lp  Lp 


L  Np  Np 


N, 


(10) 


This  method  of  describing  the  governing  equations  is  compact  and  readily 
shows  that  the  longitudinal  and  lateral-directional  equations  for  the  linearized 
case  are  uncoupled.  It  also  lends  itself  to  matrix  operations  for  determining  the 
stability  characteristics  of  the  natural  modes  for  both  longitudinal  and  lateral- 
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directional  motion.  A  general  description  of  the  natural  modes  normally 
associated  with  small  perturbation  theory  is  presented  in  Etkin  [Ref.  10].  The 
dimensional  derivatives  which  make  up  the  elements  of  and  A|at-dir  are  in 

accordance  with  NASA  convention  and  are  deHned  in  McRuer,  Ashkenas  and 
Xjraham  [Ref.  11]. 

When  the  underlying  assumptions  of  small  perturbation  theory  are  extended 
in  order  to  introduce  non*linear  terms  in  the  equations  of  motion,  the  linearized 
equations  shown  above  must  be  modified  to  account  not  only  for  the  extra  non¬ 
linear  terms,  but  also  to  drop  out  the  terms  resulting  from  the  linearization. 
Additionally,  Euler  angle  relations  are  introduced  so  that  the  roll  angle 
becomes  the  Euler  angle  0  and  the  pitch  angle  6  becomes  the  Euler  angle  6. 
The  additional  terms  to  be  added  to  transform  the  linearized  equations  into  the 
fully  coupled,  non-linear  equations  are  shown  below. 


(NL)p  =  pa +  ^(CeS,t' ■%'!>)  O') 

Iv  *  ^Z 

(NL)p  =-^1 -  qr  (12) 

*x 

(NL)0  =  (qS(])  +  iCij,)  T0  (13) 

Ix  -  ly 

(NL)r=— j-^qp  (14) 

*z 

(NL)a  =  -pP  +  u  •  %  *  V  05) 

h  *  ^x 

(NL)q  *P  +  Mct(NL)a  (16) 

(NL)0=q(C<i>-l)-rS(j)  (17) 
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and  the  new  state  vector  for  the  non-linear  equation  set  is: 


T 

Xnl  =  [P  P  O  r  a  q  ©  ]  (18) 

It  will  be  noted  that  the  aerodynamic  terms  still  retain  their  linear  form 
while  the  nonlinear  aspects  are  introduced  via  inertial  coupling.  Retaining  the 
simplified  forms  for  the  aerodynamic  terms  is  deliberate  at  this  stage  of  the 
analysis  in  order  to  illustrate  the  influence  of  inertial  coupling  upon  the  ensuing 
motion  response. 
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V.  COMPUTATIONAL  PROCEDURES 


A.  LINEAR  SYSTEMS  OF  EQUATIONS 

An  understanding  of  the  equations  of  motion  and  the  significance  of  their 
nonlinearity  is  essential  prior  to  attempting  to  determine  the  aircraft's  response 
to  an  initial  perturbation.  Had  the  equations  remained  purely  linear,  a  relatively 
simple  solution  would  have  been  available  by  calculating  the  eigenvalues  of  the 
characteristic  polynomial  and  each  associated  eigenvector.  The  general  solution 
of  a  linear  differential  equation  is  a  linear  combination  of  all  linearly 
independent  solutions.  Therefore,  the  general  solution  can  be  expressed 
explicitly  as  a  function  of  time  as  follows  (two  degree  of  freedom  system  shown 
for  clarity): 


{x(t))  *ci  (xile^^  *+C2  {X2)e^2*  (19) 

where  {x(t))  is  a  column  vector  whose  individual  components  represent  the  time 
response  of  each  degree  of  freedom,  X,i  and  \2  are  the  eigenvalues,  {xi }  and 
(x2}  are  the  associated  eigenvectors  and  ci  and  c2  are  constants  representing  the 
modal  participation  of  each  of  the  modes  in  the  response.  In  general,  the 
eigenvalues  are  complex  in  nature  and  represent  a  damped,  oscillatory  solution. 
Utilizing  the  fundamentals  of  complex  arithmetic  as  applied  to  a  response  with 
real  physical  terms,  the  expression  may  be  rewritten  as  a  combination  of  sines 
and  cosines  to  remove  the  complex  terms,  thereby  leaving  the  equation  as  a 
trigonometric  function  of  time  alone.  Once  the  time  domain  solution  is  found  in 
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this  manner,  time  history  traces  may  be  generated  using  a  short  computer 
algorithm.  Such  an  algorithm  might  iterate  time,  calculate  the  value  of  each 
component  of  (x(t)}  at  that  time,  store  the  values  in  a  data  file  and  repeat  the 
procedure.  Plotting  the  data  would  provide  the  time  history  traces.  Although 
the  procedure  outlined  here  was  intended  to  illustrate  the  ease  with  which  a 
purely  linear  system  can  be  solved,  the  technique  described  in  the  next  section 
was  used  to  solve  both  the  linear  and  nonlinear  systems  of  equations  so  that  only 
one  program  had  to  be  used  regardless  of  the  desired  form  of  the  solution.  Even 
so,  the  calculation  of  eigenvalues,  eigenvectors,  natural  frequencies  and  damping 
ratios  remained  important  and  was  used  extensively  in  the  linear  analysis  because 
the  character  of  the  dynamic  response  could  be  determined  and  visualized  solely 
by  these  parameters.  The  actual  time  history  response  provided  an  added  means 
of  visualizing  the  influence  of  these  parameters. 

B.  NONLINEAR  SYSTEMS  OF  EQUATIONS 

A  computer  program  was  used  to  numerically  integrate  both  the  linear  and 
nonlinear  equations  of  motion  to  obtain  a  time  history  response.  As  mentioned 
earlier,  the  algorithm  was  capable  of  numerically  integrating  linear  differential 
equations  as  well  as  nonlinear  equations.  Therefore,  the  same  program  was  used 
to  obtain  either  the  linear  or  the  nonlinear  response,  depending  on  the  solution 
desired. 

The  heart  of  the  program  was  a  numerical  integration  technique  based  on 
Richardson's  extrapolation;  an  explicit,  two-step  numerical  procedure  which  can 
be  readily  applied  to  first  order  differential  equations.  It  is  a  variation  of  the 
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second  order  Runge-Kutta  method.  The  basic  concept  behind  the  method  is 
presented  by  Ferziger  [Ref.  12]  and  demonstrated  by  the  following  example. 

Suppose  the  first  derivative  of  each  variable  of  interest  is  expressed 
explicitly  as  a  function  of  known  initial  conditions  at  time  t-0.  By  simple 
substitution,  the  numerical  value  of  each  derivative  can  be  calculated  and  used  to 
approximate  the  value  of  the  variables  of  interest  a  short  time  later. 
Mathematically, 


y(n+i)=yn+h  yn 


where  h  is  'he  time  increment,  n  represents  a  discretized  time  step  and  the 
superscript  (1)  denotes  that  the  value  of  y(n+l)  calculated  here  is  just  the  first 
estimate  of  the  final  value  at  time  t(n-fl)  =  tn  +  h.  The  calculation  is  then 
repeated  using  two  steps,  each  at  half  of  the  original  time  increment  as  shown 
below. 


(2)  h  . 


w  .h  ■  1/ 

y(n+l)  =  y(n+|^)+jy(n+l/2) 


The  superscript  (2)  indicates  that  the  calculated  value  is  the  second  estimate 
of  the  final  value.  The  estimates  are  then  combined  using  Richardson's 


25 


extrapolation  to  obtain  the  Hnal  expression  for  the  variable  of  interest  at  the  new 
time. 

(2)  (1) 

y(n+i)  =  2y(„+i)-y(„+i)  (23) 

This  procedure  is  repeated  for  all  of  the  variables  at  each  time  step  to  obtain 
updated  values  for  each  variable  at  the  new  time  step.  The  value  of  each  variable 
at  the  new  time  step  is  then  stored  in  a  data  file  for  subsequent  plotting.  The 
time  is  then  incremented  and  the  first  derivatives  recalculated  using  the  updated 
information  from  the  previous  time  step.  The  process  continues  for  a  length  of 
time  as  specified  by  the  user,  at  which  time  the  computer  program  terminates. 
After  an  investigation  to  determine  the  time  step  required  to  retain  sufficient 
accuracy  was  conducted,  the  time  step  was  set  at  .05  seconds. 

C.  INPUTS  TO  COMPUTER  PROGRAM 

The  program  mentioned  in  the  previous  section  required  specific  inputs 
which  served  to  identify  the  geometry,  configuration,  inertial  properties  and 
aerodynamic  characteristics  of  the  aircraft  in  its  trimmed  condition.  While  most 
of  this  information  remained  constant  regardless  of  the  angle  of  attack  studied, 
the  aerodynamic  characteristics  which  determine  the  dynamic  response  of  the 
aircraft  were  highly  dependent  upon,  among  other  things,  angle  of  attack. 
Therefore,  to  describe  the  character  of  the  aircraft  in  the  program  at  each 
different  angle  of  attack,  the  longitudinal  and  lateral-directional  "plants" 
corresponding  to  the  linear  response  at  the  desired  trim  AOA  were  used.  These 
"plants"  are  the  square  matrices  labeled  Along  and  Aiat-dir  in  equations  (06)  and 
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(10)  which  contain  the  coefficients  in  the  equations  of  motion.  These  coefHcients 
are  made  up  of  specific  combinations  of  the  stability  derivatives  which  are 
normally  determined  by  wind  tunnel  testing  and  post  flight  parameter 
identification  techniques.  The  nonlinear  terms  were  programmed  in  and  selected 
as  a  program  option  if  a  nonlinear  solution  was  chosen  by  the  user. 

In  addition  to  these  items,  initial  conditions  were  specified  to  obtain  a  non¬ 
zero  response.  The  choice  of  initial  conditions  was  important  in  obtaining  a 
coupled  response  consistent  with  the  documented  response  of  the  aircraft.  The 
choice  of  initial  conditions  will  be  discussed  later.  A  complete  listing  of  the 
computer  program  appears  in  Appendix  B. 
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VI.  ANALYSIS 


The  ultimate  goal  of  this  analysis  was  to  obtain  multiple  time  history  traces 
which  characterized  the  aircraft  wing  rock  motion,  one  trace  corresponding  to 
.each  of  the  seven  pertinent  degrees  of  freedom.  In  building  up  to  this  goal,  a 
linearized  analysis  (and,  therefore,  uncoupled  as  well)  was  used  Hrst  to  ensure 
that  the  aircraft's  response  to  known  stabilizing  or  destabilizing  conditions  would 
produce  convergent  or  divergent  behavior,  respectively.  Furthermore,  the 
linear  results  could  be  used  as  a  benchmark  to  compare  to  the  nonlinear,  coupled 
response  once  that  response  was  determined.  The  study  included  aircraft 
response  at  angles  of  attack  ranging  from  zero  to  25  degrees  so  that  the 
variations  in  the  response  due  to  the  initial  trim  condition  could  be  evaluated.  A 
summary  of  the  trim  conditions  for  flight  with  gear  up  and  flaps  at  maneuver  at 
500  ft.  appears  below  in  Table  1 . 


TABLE  1 

SUMMARY  OF  AIRCRAFT  TRIM  CONDITIONS 


fiCA 

(degrees) 

KA 

(units) 

Trim  Veiocity 
(ft/sec) 

Trim  Velocity 
(Kts) 

Mach 

number 

0 

3.3 

575 

345 

■BiM 

5 

8.8 

326 

196 

1  0 

14.3 

255 

153 

0.23 

1  5 

19.8 

228 

137 

0.20 

20 

25.3 

213 

128 

0.19 

25 

30.8 

1  94 

116 

0.17 
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The  geometric  and  inertial  properties  of  the  aircraft  described  in  the  data 
base  were  as  follows: 

WT  =  520001b. 

C.G.  =  16.2  %  MAC 
Ixx  =  51509  slug  ft^ 

Iyy  =  232773  slug  ft2 
Iz2  =  275627  slug  ft2 
Ixz=  2654  slug  ft2 

It  is  critically  important  to  recognize  that  although  an  aircraft  is  never  flown 
stick  fixed,  a  stability  analysis  conducted  in  this  manner  provides  very  useful 
information  on  the  tendencies  of  the  aircraft  to  move  about  its  trim  position  once 
disturbed.  The  flying  qualities  of  an  aircraft  are  very  much  dependent  on  these 
stability  characteristics.  In  fact,  handling  qualities  ratings  are  based  on  the 
dynamic  characteristics  of  the  linearized  modes  for  many  aircraft. 

A.  LINEAR  ANALYSIS 

To  begin  the  linear  analysis,  the  characteristic  polynomial  for  both  the 
longitudinal  and  lateral-directional  "plants"  were  solved  for  the  eigenvalues  of 
each  system.  The  natural  frequencies,  modal  damping  and  associated 
eigenvectors  were  also  calculated.  As  previously  mentioned,  these  parameters 
serve  to  identify  the  character  of  the  dynamic  response  for  each  of  the  natural 
modes.  A  summary  of  the  linearized  systems  short  period  and  dutch  roll 
characteristics  is  provide  below  in  Table  2. 
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TABLE 2 

SUMMARY  OF  LINEARIZED  DYNAMIC  CHARACTERISTICS 


Shod 

Period 

Dutch  Roll  1 

IGnTSDl 

Damping  Ratio 

Damping  Ratio 

0 

1 .8729 

0.7452 

2.5202 

0.1067 

6 

1,0607 

0.7488 

1.511 

0.0755 

10 

0.E116 

0.7352 

1.3056 

-0.0087 

15 

0.6725 

0.6642 

1.0416 

-0.1727 

20 

0.6206 

0.6716 

1.0138 

-0.3575 

25 

0.5602 

0.6766 

0.7061 

-0.4725 

A  quick  look  at  this  table  reveals  a  great  many  things  about  the  character  of 
the  aircraft's  dynamic  response  at  the  angles  of  attack  studied.  The  following 
paragraphs  provide  an  in  depth  discussion  of  the  linearized  results  for  the  F-14A. 

1.  Short  Period  Mode 

a.  Time  History  Response 

The  F'14  displays  a  stable,  heavily  damped  short  period  mode  at  all 
of  the  angles  of  attack  studied  in  this  analysis.  When  excited,  this  mode  tends  to 
produce  a  rapidly  convergent  solution  back  to  the  initial  trim  condition.  The 
short  period  response  for  zero  degrees  angle  of  attack  is  shown  below  in  Figure 
3.  The  short  period  responses  for  the  other  angles  of  attack  are  not  shown  due  to 
their  similarity  to  the  AOA=0  response. 
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Figure  3.  Short  Period  Response  for  AOA  *  0  degrees 


The  initial  conditions  used  to  generate  the  short  period  response 
corresponded  to  a  one*tenth  scaling  of  the  normalized  complex  eigenvector  at 
t=0.  The  normalized  eigenvector  was  scaled  in  this  muiner  to  produce  an  initial 
perturbation  quantity  consistent  with  the  magnitudes  expected  in  flight.  This 
technique  took  into  account  the  phase  relationship  between  each  component  of  the 
eigenvector  and  produced  the  true  physical  relationship  between  each  component 
at  any  time. 

b.  Root  Locus 

It  is  interesting  to  note  the  changing  position  of  the  characteristic 
roots  (eigenvalues)  of  the  longitudinal  system  as  angle  of  attack  changes.  The 
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migration  of  the  roots  has  a  large  impact  on  the  character  of  the  dynamic 
response  in  that  the  natural  frequency  and  damping  ratio  are  determined  by  the 
location  of  the  roots  in  the  complex  plane.  Shown  below  in  Figure  4  is  a  root 
locus  showing  the  migration  of  the  short  period  roots. 


2.  Dutch  Roll  Mode 

The  F-I4  displays  an  unstable  dutch  roll  mode  at  high  angle  of  attack 
due  in  part  to  degraded  directional  stability  (Cnp  going  positive  to  negative) 
which  appears  at  around  15-20  degrees  AOA,  depending  on  aircraft  speed  and 
external  stores  loading.  The  Out  of  Control  Flight  Training  Guide  {Ref.  13] 
illustrates  this  degradation  with  speed  and  external  loading  as  shown  in  Figures  5 
and  6.  The  variation  of  the  directional  stability  parameter  Cnp  calculated  from 
the  data  base  is  shown  in  Figure  7.  The  increase  in  dihedral  effect  (C|p  getting 
more  negative)  and  the  decrease  in  roll  damping  (Cjp  getting  less  negative)  as 
angle  of  attack  increases  also  contribute  significantly  to  the  dutch  roll  instability. 
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DIRECTIONAL  STABIUTY 
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Figure  S.  Influence  of  Speed  on  Directional  Stability 
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F-14  EXHIBITS  NEGATIVE  DIRECTIONAL  STABILITY  AT  HIGH  AOA 

AT  LOW  AOA,  ADVERSE  EFFECT  CAN  BE  COUNTERED  BY  DIHEDRAL  EFFECT 

AT  HIGH  AOA.  POSITIVE  PILOT  ACTION  REQUIRED  TO  RECOVER 


Alpha  (degrees) 


Figure  7.  Variation  of  Cnp  with  Angle  of  Attack  from  Database 

When  this  unstable  mode  is  excited,  the  aircraft  tends  to  diverge  unless 
pilot  inputs  are  made  to  control  the  motion. 
a.  Time  History  Response 

The  linear  analysis  reveals  that  the  dutch  roll  mode  is  lightly 
damped  at  low  angles  of  attack.  This  results  in  a  slowly  convergent  oscillation 
for  both  zero  and  five  degrees  angle  of  attack.  As  angle  of  attack  increases  past 
ten  degrees,  the  dutch  roll  mode  goes  unstable.  At  this  point  it  is  only  slightly 
unstable,  resulting  in  a  very  slow  divergent  oscillation.  As  angle  of  attack 
continues  to  increase,  the  dutch  roll  mode  becomes  extremely  unstable,  resulting 
in  a  rapid  divergence.  Figures  8,  9  and  10  show  the  dutch  roll  time  history 
responses  for  AOA=0,  10  and  15  degrees,  respectively. 
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Dutch  Roll  Response  at  AOA  =  10  degrees 


The  initial  conditions  for  all  of  the  dutch  roll  responses  used  a  one- 
tenth  scaling  of  the  normalized  complex  eigenvector,  just  as  in  the  short  period 
case. 

b.  Root  locus 

An  explanation  for  the  behavior  of  the  dutch  roll  mode  is  clear 
when  an  examination  of  the  placement  of  the  roots  in  the  complex  plane  is  made. 
A  root  locus  is  shown  below  in  Figure  11  for  the  lateral-directional  system's 
dutch  roll  roots.  Notice  the  migration  of  the  characteristic  roots  through  the 
imaginary  axis  as  the  angle  of  attack  approaches  ten  degrees.  Once  the  roots 
migrate  to  the  right  hand  half  plane,  the  response  changes  from  a  convergent  to  a 
divergent  oscillation. 

4 

2 

o 
<u 
(0 

-2 

~4 

-0.4  -0.3  -0.2  -0.1  0  0.1  0.2  0.3  0.4 

Figure  1 1 .  Dutch  Roll  Root  Locus 

As  previously  stated,  the  instability  of  the  dutch  roll  mode  at  high 
angle  of  attack  results  primarily  from  changes  in  Cnp,  Cip  and  CIp.  The  root 
locus  plots  in  Figures  12,  13  and  14  show  the  sensitivity  of  the  roll,  spiral  and 
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dutch  roll  roots  to  individual  changes  in  each  of  these  parameters.  The  direction 
of  root  migration  shown  corresponds  to  variation  of  the  stability  parameter  from 
the  value  at  AOA*0  to  AOA=25.  All  other  parameters  correspond  to  AOAa25 
values.  When  the  changes  in  CnP,  Cip  and  CIp  occur  simultaneously  as  AOA 
increases,  the  individual  effects  combine  to  produce  a  rapid  onset  of  dutch  roll 
instability  as  shown  previously  in  Hgure  II.  It  was  noted  that  had  the  values  of 
CnP,  Cip  and  CIp  corresponding  to  zero  AOA  been  held  constant  as  AOA 
increased,  the  dutch  roll  roots  would  have  been  stable. 


Real  (rad/sec) 

Figure  12.  Root  Locus  for  Cpp  Variation 
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Figure  13.  Root  Locus  for  Cjp  Variation 
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Figure  14.  Root  Locus  for  Qp  Variation 


B.  NON-LINEAR  ANALYSIS 

Now  that  an  understanding  of  the  effect  of  angle  of  attack  on  the  linearized 
system  dynamics  has  been  gained,  a  comprehensive  study  of  the  nonlinearized 
system  dynamics  is  in  order.  This  begins  by  selecting  the  coupling  option  in  the 
computer  program,  which  has  the  effect  of  including  the  nonlinear  terms  in  the 
equations  of  motion  before  numerical  integration  takes  plac^^.  Additionally, 
initial  conditions  are  chosen  to  demonstrate  the  effect  of  coupling  from  the 
lateral-directional  system  to  the  longitudinal  system.  Specifically,  the  lateral- 
directional  initial  conditions  are  the  now  familiar  scaled  eigenvector  components, 
while  the  longitudinal  initial  conditions  are  set  to  zero. 

At  the  low  angles  of  attack  where  both  linear  modes  are  convergent,  stable 
oscillations,  one  would  expect  that  coupling  of  the  two  systems  of  equations 
would  also  produce  a  convergent,  stable  solution.  This  was  exactly  the  case  and 
is  clearly  evident  from  examination  of  Figure  15,  the  nonlinear  response  at  zero 
degrees  angle  of  attack.  The  changes  in  the  lateral-directional  parameters  from 
the  linear  dutch  roll  response  were  imperceptible  at  this  angle  of  attack.  The 
effect  of  coupling  on  the  longitudinal  parameters  is  also  evident  from  the  plots, 
although  the  maximum  amplitude  of  the  longitudinal  perturbations  is  on  the 
order  of  a  tenth  of  a  degree.  Although  the  longitudinal  response  shown  here  is 
imperceptible  from  a  pilot's  point  of  view,  there  is  some  significance  to  one 
feature  which  continues  to  appear  at  higher  angles  of  attack.  This  is  the 
frequency  relationship  between  the  roll  angle  response  and  the  angle  of  attack 
response.  In  every  case,  the  angle  of  attack  perturbations  occurred  at  twice  the 
frequency  of  the  roll  response.  This  is  somewhat  intuitive,  as  angle  of  attack 
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Figure  IS.  Coupled  Response  at  AOA  =  0  degrees 
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perturbations  can  be  expected  to  occur  with  either  right  or  left  wing  movement. 
This  tendency  has  been  observed  in  the  actual  flight  test  time  history  traces 
discussed  earlier. 

A  similar  response  was  obtained  for  five  degrees  AOA.  Once  again,  the 
nonlinear  coupling  produced  very  little  change  from  the  linear  dutch  roll 
response  and  extremely  small  longitudinal  perturbations.  At  ten  degrees  AOA, 
however,  the  nonlinear  response  failed  to  converge.  Instead,  a  stable,  constant 
magnitude  oscillation  developed  in  the  lateral-directional  parameters  which  is 
characteristic  of  a  mild  wing  rock  motion..  This  response  is  shown  in  Figure  16. 
The  maximum  amplitudes  of  the  longitudinal  parameters  due  to  the  coupling  are 
increasing,  but  are  still  insignificant  in  comparison  to  the  overall  response.  Note 
that  the  mean  value  of  the  longitudinal  perturbations  in  angle  of  attack  and  pitch 
rate  are  displaced  from  their  respective  equilibrium  positions,  while  the  pitch 
angle  begins  to  fall  off  in  an  oscillatory  manner.  The  presence  of  limit  cycle 
behavior  is  clearly  shown  in  this  figure. 
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Figure  16.  Coupled  Response  at  AO  A  =  10  degrees 


Beyond  ten  degrees  angle  of  attack,  the  coupling  continued  to  produce 
oscillations  in  the  lateral-directional  paranteters  which  were  indicative  of  wing 
rock,  while  more  dramatic  changes  began  to  appear  in  the  longitudinal 
parameters.  Responses  for  AOA:=15,  20  and  25  degrees  are  shown  on  the 
following  pages  in  Figures  17  through  19. 

The  results  of  this  portion  of  the  analysis  further  demonstrated  the 
development  of  wing  rock  limit  cycles  as  a  result  of  nonlinear  coupling  at  high 
angle  of  attack.  Although  the  maximum  amplitudes  for  the  roll  angle  and 
sideslip  were  somewhat  higher  than  those  given  in  NATOPS  for  the  maneuver 
flap  conHguration,  the  magnitudes  were  not  unreasonable  considering  the 
approximations  made  to  obtain  the  results. 

A  few  other  important  observations  can  be  made  by  referring  back  to  the 
plots  of  sideslip,  roll  angle,  angle  of  attack  and  pitch  angle  for  AOA=20,  which 
have  been  enlarged  in  Figures  20a  and  20b  to  show  greater  detail.  The 
nonlinearity  of  the  sideslip  response  near  the  maximum  amplitude  limits  of  the 
oscillation  can  be  seen  in  Figure  20a.  Instead  of  obtaining  the  smooth  sinusoidal 
response  normally  associated  with  an  underdamped  system,  this  response 
illustrates  that  the  coupling  distorts  the  oscillation.  Second,  it  is  noted  that  as  the 
perturbations  in  sideslip  and  roll  angle  begin,  the  frequency  of  the  oscillations  is 
very  close  to  the  dutch  roll  frequency.  As  the  coupling  takes  effect  and  the  limit 
cycle  is  established,  the  frequency  of  oscillation  decreases  by  approximately  30 
percent.  This  trait  was  observed  at  all  AOA's  where  a  limit  cycle  was  established 
and  can  be  attributed  to  the  nonlinear  interaction  of  the  longitudinal  and  lateral- 
directional  modes.  For  the  aircraft  modeled  in  Ref.  4  which  displayed  a  very 
mild  dutch  roll  instability,  this  phenomenon  was  not  observed,  In  that  study,  the 
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frequencies  of  the  natural  modes  were  artiHcially  controlled  to  maintain  a 
harmonic  relationship  between  the  short  period  and  the  dutch  roll  response. 
Here,  no  such  relationship  exists.  Further  study  may  uncover  the  dependency  of 
the  limit  cycle  frequency  on  the  characteristic  frequencies  of  the  natural  modes. 

Finally,  it  was  noted  that  as  the  limit  cycle  was  established,  the  mean  value  of 
the  angle  of  attack  perturbation  increased  to  a  value  which  was  slightly  higher 
than  its  equilibrium  position  and  the  pitch  angle  began  to  drop.  This  tendency 
was  also  observed  in  the  flight  test  results.  Although  the  reasons  for  this  are 
complex  as  a  result  of  multiple  dependency  between  the  perturbation  quantities  in 
the  fully  coupled  nonlinear  equations,  one  aspect  of  the  wing  rock  motion  may 
provide  a  clue  as  to  its  origin.  Consider  the  lift  equilibrium  during  the  wing 
rock  motion.  The  vertical  component  of  lift  at  any  time  is  proportional  to  the 
cosine  of  the  roll  angle,  which  changes  continuously  with  time.  Therefore,  the 
average  value  of  lift  during  one  limit  cycle  oscillation  is  dependent  upon  the 
maximum  amplitude  of  the  wing  rock  and  will  always  be  lower  than  the 
equilibrium  value  of  lift.  As  a  result  of  the  drop  in  the  average  lift  per  cycle,  the 
aircraft  starts  to  travel  on  a  curvelinear  flight  path.  This  response  shows  up  as  a 
loss  in  altitude.  Additionally,  the  dependence  of  the  nonlinear  terms  in  the  angle 
of  attack  and  pitch  angle  equations  ((15)  and  (17))  upon  the  average  value  of  lift 
is  readily  apparent  by  the  presence  of  the  cos(^)  term.  This  dependency  may  be 
the  dominant  feature  w'hich  causes  this  type  of  response. 
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VII.  RESULTS 


The  convergent  nature  of  the  linearized  short  period  response  and  the 
divergent  nature  of  the  linearized  dutch  roll  response  for  the  same  flight 
condition  indicated  that  the  numerical  procedure  used  produced  time  history 
traces  consistent  >vith  the  behavior  of  stable  and  unstable  linear  systems, 
respectively. 

When  the  equations  were  modified  to  include  the  nonlinear  terms,  the 
responses  for  the  low  angles  of  attack  (AOA=0  and  5  degrees)  did  not  change 
appreciably,  although  the  influence  of  the  coupling  was  apparent  in  small 
longitudinal  perturbations.  At  these  angles  of  attack,  both  the  longitudinal  and 
lateral-directional  characteristic  roots  were  stable,  resulting  in  an  asymptotically 
stable  solution  in  the  sense  of  Liapunov  [Ref.  14].  It  has  been  demonstrated  that 
although  coupling  between  the  two  sets  of  equations  occurred,  the  overall  result 
failed  to  produce  a  limit  cycle  response. 

When  the  nonlinear  equations  were  analyzed  at  the  higher  angles  of  attack, 
however,  the  presence  of  unstable  dutch  roll  eigenvalues  destroyed  the 
asymptotic  stability.  Even  so,  the  coupling  between  the  longitudinal  and  lateral- 
directional  equations  still  yielded  a  stable  system  response  in  the  form  of  a  wing 
rock  limit  cycle  oscillation. 
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VIII.  CONCLUSIONS 

A  numerical  an'Uysis  of  the  nonlinear  equations  of  motion  has  been 
conducted  to  investigate  the  contribution  of  inertial  coupling  to  the  development 
of  wing  rock  in  the  F*14  aircraft.  Actual  wind  tunnel  data  was  used  to  develop 
all  of  the  stability  parameters  for  the  analysis.  Although  a  number  of 
simplifying  assumptions  were  made,  the  analysis  indicated  that  inertial  coupling 
of  a  stable  short  period  mode  and  an  unstable  dutch  roll  mode  can  result  in  a 
response  very  much  like  that  encountered  in  the  aircraft,  especially  for  the 
lateral-directional  parameters.  The  trends  displayed  in  the  longitudinal 
parameters  as  a  result  of  the  coupling  were  consistent  with  flight  tests;  however, 
the  magnitude  of  the  excursions  from  the  trim  position  for  these  parameters  was 
greater  than  expected.  Although  these  deviations  appeared  to  be  beyond  reason, 
the  results  obtained  here  should  not  be  discounted  on  this  basis.  The  intention 
was  to  conduct  a  preliminary  investigation  into  the  mechanics  of  wing  rock  for 
the  F- 14  in  hopes  of  uncovering  a  relatively  simple  explanation  for  the  aircraft's 
behavior.  In  the  pursuit  of  tliis  goal,  it  has  been  demonstrated  that  a  stable  short 
period  mode  can  feed  damping  energy  into  an  unstable  dutch  roll  mode  to 
produce  a  bounded  wing  rock  type  oscillation.  Unquestionably  there  is  a  great 
deal  of  further  research  which  must  be  done  to  fully  understand  the  motion, 
including  more  detailed  analyses  which  disregard  some  of  the  simplifications 
made  in  this  analysis  and  which  promise  to  yield  results  which  are  more 
consistent  with  the  aircraft's  actual  response.  The  following  section  highlights  a 
few  of  these  techniques. 


IX.  RECOMMENDATIONS  FOR  FURTHER  RESEARCH 


The  analysis  presented  here  is  only  a  start.  It  serves  to  illustrate  that 
numerical  techniques  can  be  used  to  solve  the  nonlinear  equations  of  motion  and 
to  predict  the  response  of  aircraft  subject  to  specified  initial  conditions.  There 
exists  a  wide  variety  of  related  research  topics  that  can  further  our  understanding 
of  nonlinear  flight  mechanics.  Some  of  these  topics  are  addressed  in  the 
paragraphs  which  follow. 

A.  EIGHT  DEGREE  OF  FREEDOM  ANALYSIS 

Perhaps  the  least  difficult  of  all  suggestions  to  follow  for  continued  research 
in  this  field  would  be  the  extension  of  the  analysis  to  eight  degrees  of  freedom. 
The  eighth  parameter  of  interest,  pcnurbation  velocity,  was  held  constant  during 
this  analysis  based  on  an  approximation  of  constant  velocity  during  the  wing  rock 
motion.  To  bring  the  analysis  to  eight  degrees  of  freedom,  the  X-force  equation 
would  be  introduced  in  its  full  non-linear  form,  the  pertinent  stability  parameters 
obtained  from  further  analysis  of  the  available  data  base  and  the  computer  code 
modified.  The  influence  of  the  additional  degree  of  freedom  on  both  the 
linearized  analysis  of  the  longitudinal  parameters  and  the  full  non-linear  analysis 
would  wanant  a  full  investigation  of  the  changes  in  the  dynamic  response, 
especially  given  the  tendency  of  the  aircraft  to  end  up  in  a  nose  down  attitude 
where  changes  in  velocity  are  sure  to  occur. 
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B.  TIME  DEPENDENT  STABILITY  PARAMETER  ANALYSIS 

A  much  more  difficult  undertaking  would  involve  the  incorporation  of  time 
dependent  stability  parameters  into  the  analysis.  As  stated  previously,  the  data 
base  contains  an  extensive  list  of  stability  parameters  as  a  function  of  multiple 
variables;  however,  only  derivative  values  corresponding  to  the  steady,  level 
flight  trim  condition  were  used  for  this  investigation.  Recall  that  as  the  aircraft 
moves,  the  pitch,  roll  and  yaw  rates  developed  contribute  to  the  stability  of  the 
aircraft  and  may  play  a  significant  role  in  the  overall  response  of  the  aircraft. 
Furthermore,  changes  in  AOA  and  sideslip  angle  as  the  aircraft  moves  about  its 
equilibrium  position  also  affect  the  stability  parameters.  Lastly,  although  this 
analysis  was  conducted  controls  fixed,  ..  :  actual  aircraft,  control  surface 

position  also  influences  the  stability  c.  aircraft.  Therefore,  although  this 
analysis  proceeded  under  the  approximation  of  constant  stability  characteristics 
as  part  of  a  planned  buildup  to  a  more  detailed  analysis,  it  nevertheless  neglects 
the  changing  dynamic  behavior  of  the  aircraft  as  it  translates  and  rotates.  Ease 
of  implementation  for  the  first  attempt  at  numerical  creation  of  the  wing  rock 
motion  from  wind  tunnel  data  was  also  a  factor  in  choosing  this  approximation. 
Although  this  simplification  may  be  appropriate  for  some  aircraft  operating  at 
relatively  low  AOA,  it  may  be  the  reason  for  obtaining  such  large  perturbations 
in  the  longitudinal  parameters  at  high  AOA. 

Ii.  order  to  implement  such  a  change  in  the  analysis,  one  of  several  options 
must  be  exercised  to  allow  for  variable  stability  parameters.  One  option  might 
involve  the  development  of  subroutines  which  would  allow  the  computer  to 
conduct  a  table  look  up  of  dimensionless  stability  parameters  and  subsequent 
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conversion  into  dimensional  form.  This  would  require  the  availability  of  the 
entire  data  base  to  the  code,  as  well  as  restructuring  of  the  code  to  accommodate 
iflie  table  look  up  feature. 

Another  option  might  involve  the  use  of  multi-variable  curve  fitting  to 
obtain  very  accurate  approximations  of  the  dimensionless  stability  parameters  as 
the  aircraft  moves. 

C.  OPTIMIZATION  OF  NUMERICAL  SOLUTION 

Still  another  possibility  for  future  research  in  this  area  involves  numerical 
optimization.  Numerical  techniques  commonly  used  to  obtain  solutions  for 
systems  of  differential  equations  vary  widely  in  their  accuracy  and  suitability  for 
a  given  problem.  The  "exactness*'  of  a  numerical  solution  is  very  much 
dependent  upon  the  methods  used  to  numerically  approximate  the  equations,  the 
desired  accuracy  of  the  solution  and  the  acceptable  numerical  cost  in  detennining 
that  solution.  Similarly,  the  stability  of  a  numerical  procedure  is  also  quite 
important  in  that  it  determines  the  acceptable  variation  of  parameters,  such  as  the 
time  increment,  v/hich  force  a  convergent  solution.  Care  must  be  taken  when 
conducting  numerical  analysis  to  account  for  the  stability  characteristics  and 
desired  accuracy  of  the  solution.  Therefore,  a  study  attempting  to  optimize  the 
numerical  technique  used  to  obtain  the  time  history  data  points  may  be  in  order. 
For  example,  a  number  of  different  finite  difference  schemes  could  be  used  to 
obtain  and  compare  solutions,  computation  time,  overall  efficiency,  etc. 

D.  INCORPORATION  OF  ACTUAL  FLIGHT  TEST  RESULTS 

Once  a  numerical  technique  is  developed  to  more  accurately  and  efficiently 

model  all  aspects  of  the  wing  rock  motion,  the  incorporation  of  additional  flight 
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test  results  from  an  instrumented  aircraft  could  be  used  to  verify  results  obtained 
from  the  study.  By  utilizing  the  same  initial  conditions,  aircraft  configuration, 
geometry  and  inertia  characteristics  as  inputs  into  a  similar  analysis,  it  would  be 
possible  to  try  to  numerically  reconstruct  the  true  response  of  the  aircraft. 

E.  NUMERICAL  ANALYSIS  OF  F/A-18  WING  ROCK 

In  that  a  complete  F/A~18  flight  simulation  program  and  data  base  are 
available  for  research  use,  a  similar  analysis  conducted  on  that  aircraft  may 
provide  some  clues  as  to  the  mechanics  behind  F/A-18  wing  rock.  It  is  known 
that  the  aircraft  displays  a  complex  wing  rock  motion  at  very  high  AOA  which 
cannot  be  damped  out  by  pilot  input,  flight  control  computer  or  any  combination 
of  the  two.  The  data  base  can  be  obtained  from  the  NASA  Ames  Research 
Center/Dryden  Flight  Research  Facility  located  at  Edwards  Ara,  CA. 
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APPENDIX  A.-  DATA  BASE  MANIPULATION 

j 

The  stability  derivatives  needed  to  conduct  this  analysis  were  extracted  from 
tabular  wind  tunnel  data  provided  by  the  Grumman  Corporation  [Ref.  9].  This 
-  appendix  is  provided  to  give  the  reader  some  insight  as  to  how  the  tabular  data 
was  used  to  extract  the  necessary  information  for  the  study. 

The  first  things  needed  to  conduct  the  analysis  at  any  particular  angle  of 
attack  were  the  velocity,  thrust  required  and  control  surface  positions  for  steady, 
level  flight  (SLF)  at  that  AOA.  Assumptions  of  perfect  lateral  symmetry  about 
the  X-Z  plane  and  the  absence  of  gyroscopic  effects  led  to  the  assignment  of  zero 
deflection  for  any  lateral  or  directional  control  surfaces.  The  trim  velocity  and 
stabilator  position  were  found  by  simultaneously  solving  the  following  two 
equations  for  alpha  and  Ss  (stabilator  position)  in  an  iterative  process: 

*^cg°*^0'^‘^lhrusl*‘^a“  **^85**  (23) 

CL  =  CL(,  +  CL,^^^^,  +  CL^a  +  CLjfe  (24) 

Specifically,  the  desired  trim  AOA  was  designated  and  a  guess  at  the 
appropriate  trim  velocity  was  made.  Cmcg  zero  in  SLF  and  Cl  can  be 
calculated  easily  by  using  the  assumed  trim  velocity  and  the  weight  of  the 
aircraft.  Values  for  CmO,  Cma.  CtOf  and  Clo  were  obtained  by  constructing 
crossplots  of  the  data  given  in  the  tables  to  get  the  appropriate  slope.  Local 
derivatives  were  used  in  all  cases  where  the  data  was  nonlinear.  CmSs  was  taken 
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directly  from  the  tables  and  CtSs^as  obtained  by  doing  a  coordinate 
transformation  through  the  angle  of  attack  of  Cx5g  and  CzSs  data  in  the  table.  In 
order  to  find  CLihrust  and  Cmthnisi,  the  thrust  required  for  SLF  had  to  be 
determined  by  finding  the  drag  in  SLF  at  the  desired  AOA.  This  was 
accomplished  by  doing  another  coordinate  transformation  on  the  C^p.  Czp>  Qss 
and  CzS,  information  in  the  table.  A  plot  of  Cd  vs  alpha  was  constructed  from 
this  data  for  each  different  stabilator  setting  given  in  the  tables.  From  this  plot, 
an  estimate  of  the  drag  could  be  determined.  This  plot  is  shown  below  in  Figure 
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Figure  21.  Cd  vs  Alpha 
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With  an  estimate  of  the  drag,  thrust  required  is  immediately  known.  This  is 
converted  into  a  pitching  moment  contribution  and  a  lift  contribution,  which  are 
then  converted  into  their  final  coefficient  form.  With  all  of  the  parameters  for 
these  two  equations  known,  the  equations  are  solved  for  alpha  and  8s.  The  initial 
assumed  trim  velocity  is  iterated  until  a  match  is  made  between  the  calculated 
AOA  and  the  desired  AOA.  When  a  match  is  made,  the  stabilator  setting 
calculated  is  used  to  update  the  determination  of  drag.  The  equations  are  solved 
again  with  this  new  drag  figure  and  after  only  a  few  iterations,  the  calculation 
converges  on  the  final  result  for  stabilator  setting. 

A  convenient  cross  check  for  this  procedure  is  available  by  building  up  a  plot 
of  Cmcg  vs.  alpha  by  combining  Cmo.  CmSs  and  Cmihrusi  information.  This 
procedure  results  in  the  plot  shown  below  in  Figure  22.  This  plot  can  be  used  to 
determine  the  stabilator  setting  required  to  maintain  Cmcg  =  0  O-c.,  SLF)  at  any 
desired  AOA. 

The  determination  of  stabilator  setting  is  important  because  some  of  the  data 
in  the  tables  is  dependent  on  this  parameter.  Once  this  setting  is  known,  all  of 
the  stability  derivatives  can  be  found.  The  following  steps  contain  a  brief 
description  of  the  remaining  procedures  to  obtain  these  derivatives. 

1)  Cyp,  Cnp,  C|p  and  Cza  are  found  by  constructing  crossplots  of  the  data 
given  in  the  tables  and  taking  local  derivatives. 

2)  Cyp,  Cyr,  Cip,  Cif,  Cnp,  Cnr  and  Cmq  are  taken  directly  from  the  tables. 

3)  Clq  is  found  by  doing  a  coordinate  transformation  of  C^q  and  Czq. 
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All  of  these  dimensionless  coefficients  are  then  converted  into  dimensional 
form  using  standard  conversion  formulas.  At  this  point,  the  derivatives  are 
ready  to  be  used  in  the  analysis. 


Cmcg  vs .  Alpha 
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Figure  22,  C^cg  vs  Alpha 
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APPENDIX  B.-  COMPUTER  PROGRAM 


100  REM  F14  WING  ROCK  PROGRAM 

110  REM  DATA  FOR  THE  FI  4  AT  AOA-20  DEGREES,  Vtrim-213  fp8 

120  REM  A( . )  -  PLANT  MATRIX 

130  REM  XK( )  .  PERTURBATION  QUANTITY  (PQ) 

140  REM  XK1  ( )  -  ESTIMATION  OF  PQ  AT  NEXT  FULL  TIME  STEP  BASED 

1 50  REM  ON  DERIVATIVE  AT  CURRENT  TIME  STEP 

160  REM  XK15( )  -  ESTIMATION  OF  PQ  AT  HALF  TIME  STEP  BASED  ON 

170  REM  DERIVATIVE  AT  CURRENT  TIME  STEP 

180  REM  XK2( )  -  ESTIMATION  OF  PQ  AT  FULL  TIME  STEP  BASED  ON 

1 90  REM  DERIVATIVE  AT  HALF  TIME  STEP 

200  REM  XSUM( )  -  DERIVATIVE  AT  HALF  OR  FULL  TIME  STEP 

210  REM  NL( )  -  NONLINEAR  DERIVATIVE  TERMS 

220  REM  XKMAX( )  -  MAX  VALUE  OF  PQ 

300  DIM  A(10.10),XK(10),XK1(10).XK15(10),XK2(10) 

310  DIM  XSUM(10).NL{10),XKMAX(10) 

320  ALPHA-20 
330  VTHIM-213 
340  PI-3,1415927# 

350  U-VTRIM‘COS(ALPHA*PI/180) 

360  REM  DEFINE  INERTIAL  PROPERTIES 

370  AIX-51 509l:AIY-232773l:AIZ-275627l:AIX2-2654 

380  REM  DEFINE  TERMS  IN  EQNS  OF  MOTION 

390  AIQRL-(AIY-AIZ)/AIX  :  AlPQN-  (AIX.AIY)/AIZ 

400  AIRPM-(AIZ-AIX)/AIY 

410  AIPBM  -  -,16  :  AZADOT-0  :THET0  -ALPHA*PI/1 80 
420  REM  AIPBM  is  Malphadot,  AZADOT-Zalphadot 
430  GU  -  32.1 7/U 
440  DT-,05  :  DT2-.5*DT 

500  REM  DATA  FOR  PlANTS.  ROWS  1-4  LAT-DIR,  ROWS  5-7  LONG 
510  DATA  -.0491, .0035, .151 1,-1.0007,0.0.0 
620  DAT/  -8.338,-, 6290,0, ,6877.0.0.0 
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530  DATA  0,1 ,0,0, 0,0,0 
540  DATA  -.0515,-.0692,0,-.1186,0,0,0 
550  DATA  0,0,0,0,-.2671,.9659,-.0550 
-580  DATA  0,0,0,0,-.2285,-.5278,.0088 
570  DATA  0,0,0,0,0,1,0 
580  REM  READ  IN  PLANT  DATA  TO  MATRIX  A( ) 

590  FOR  I-1  TO  7 
s600  FOR  J-1  TO  7 
610  READ  A(I,J) 

620  NEXT  J 
630  NEXT  I 

700  REM  INITIALIZE  MATRICES  FOR  COMPUTATION 
710  FOR  1-1  TO  7 

720  XK(l)-0  :  XK1(l)-0  :  XK15(l)-0  :  XK2(l)-0 
730  NEXT  I 

800  REM  DEFINE  1/10th  SCALE  DR  EIGENVECTOR  IC’S 
810  XK(1)-  -.01434:  XK(2).0:  XK(3)-.09212  :  XK(4)-.0055 
820  REM  SET  SP  ICS  TO  ZERO 
830  XK(5)-  0:  XK{6)«0:  XK(7)-0 

900  REM  INPUT  "USE  INERTIAL  COUPLING  (Y/N)";A$ 

910  REM  INPUT  "CREATE  DATA  FILES  (Y/N)":B$ 

920  A$-"Y" ;  B$-"Y" 

930  IF  B$-"Y"  THEN  GOSUB  2000 
940  REM  PRINT  INITIAL  VALUES  OF  PQ’S 
950  GOSUB  2100 

1000  REM  PERFORM  NUMERICAL  INTEGRATION 
1010  T-0 

1020  FOR  K-1  TO  400 
1030  T-K‘DT 

1 1 00  REM  CALCULATE  DERIVATIVES  AT  CURRENT  TIME  STEP 
1 1 00  REM  USING  PQ'S  AT  CURRENT  TIME  STEP 
1120  FOR  1-1  TO  7 
1130  XSUM(l)-0 


65 


1140  F0RJ-1T0  7 

1150  XSUM(I)-XSUM(I)+A(I,J)*XK(J) 

1160  NEXTJ 

1 200  REM  CALCULATE  ESTIMATED  PQ’S  AT  HALF  AND  FULL  TIME 
1 21 0  REM  STEP  BASED  ON  DERIVATIVES  AT  CURRENT  TIME  STEP 
1220  XK15(I)-XK(I)+DT2*XSUM(I) 

1230  XK1(I)-XK(I)+DT*XSUM(I; 

1240  NEXT  I 

1 250  REM  IF  NECESSARY.  CALCULATE  NONLINEAR  TERMS 
1260  IF  A$-"N"  THEN  GOTO  1420 
1270  GOSUB  2200 

1 280  REM  CORRECT  XKl 5( )  &  XK1  ( )  FOR  NONLINEAR  TERMS 

1290  FOR  1-1  TO  7 

1300  XK15(I)-XK15(I)+DT2*NL{I) 

1310  XK1(I)-XK1(I)+DT*NL(I) 

1320  NEXT  I 

1400  REM  CALCUALTE  DERIVATIVE  AT  HALF  TIME  STEP 

1410  REM  USING  ESTIMATED  PQ'S  AT  HALF  TIME  STEP 

1420  FOR  1-1  TO  7 

1430  XSUM(l)-0 

1440  FOR  J«1  TO  7 

1 450  XSUM(I)-XSUM(I)+A(I.J)*XK1 5(J) 

1460  NEXT  J 

1 500  REM  CALCULATE  ESTIMATED  PQ’S  AT  FULL  TIME  STEP 
1 51 0  REM  BASED  ON  DERIVATIVE  AT  HALF  TIME  STEP 
1520  XK2(I)»XK15(I)+DT2*XSUM(I) 

1530  NEXT  I 

1 540  REM  IF  NECESSARY.  CALCULATE  NONLINEAR  TERMS 
1550  IF  A$-"N*'  THEN  GOTO  1720 
1560  GOSUB  2400 

1 570  REM  CORRECT  XK2( )  FOR  NONLINEAR  TERMS 

1580  FOR  I-  1  TO  7 

1590  XK2(I)-XK2(I)+DT2*NL(I) 

1600  NEXT  I 
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1 700  REM  APPLY  RICHARDSON’S  EXTRAPOLATION  TO  GET 
1710  REM  FINAL  VALUE  AT  NEW  TIME  STEP 
1720  FOR  1-1  TO  7 
1730  XK(I)-2*XK2(I)-XK1(I) 

1 800  REM  FIND  MAX  VALUE  OF  EACH  PQ 

1810  IF  ABS(XK(l))>XKmax(l)  THEN  XKmax(l)-ABS(XK(l)) 

1900  NEXT  I 

1910  REM  PRINT  RESULTS 

1920  IF  (INT(K/4).(K/4))  THEN  GOSUB  2100 

1930  NEXT  K 

1940  IF  B$-"Y"  THEN  GOSUB  2700 
1950  END 

2000  REM  S/R  FOR  OPENING  DATA  FILES 
2010  OPEN  "data"  FOR  OUTPUT  AS  #1 
2020  OPEN  "maxdata"  FOR  OUTPUT  AS  #2 
2030  RETURN 

2100  REM  S/R  PRINTS  TIME  AND  STATE  VECTOR  ON  SCREEN 
2110  REM  AND  STORES  DATA  IN  FILE  FOR  SUBSEQUENT  PLOTTING 
2120  PRINT  USING  "###.####":T;XK(1):XK(2):XK(3);XK(4), 

I  XK(5).XK(6).XK(7) 

2130  IF  B$-"Y"  THEN  GOSUB  2600 
2140  RETURN 

2200  REM  S/R  FOR  FIRST  USE  OF  NON  LINEAR  TERMS 
2210  THETA  -THET0+XK(7)  :  CTHETA-COS(THETA)  : 

!  TTHETA-TAN(THETA) 

2220  PHI-XK(3)  :  SPHI-SIN(PHI)  :  CPHI-COS(PHI) 

2230  NL(1)-XK(2)‘XK(5)+GU*CTHETA*SPHI*A(1.3)‘XK(3) 

2240  NL(2)-AIQRL"XK(6)*XK(4) 

2250  NL(3)-(XK(6)*SPHI+XK(4)*CPHI)*TTHETA 
2260  NL(4)-AIPQN‘XK(2)*XK(6) 

2270  NL(5)-XK(2)*XK{1)*U/(U-AZADOT)+32.17/(U-AZADOT) 

!  *(CTHETA*CPHI-COS(THETO)+SIN(THETO)*XK(7)) 

2280  NL(6)-AIRPM*XK(4)‘XK(2)+AIPBM*NL(5) 
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2290  NL(7)-XK(6r(CPHI-1  )-XK(4)‘SPHI 
2300  RETURN 

2400  REM  S/R  FOR  SECOND  USE  OF  NON  LINEAR  TERMS 
2410  THETA  -THET0+XK15(7) :  CTHETA-COS(THETA) : 

!  TTHETA-TAN(THETA) 

2420  PHI-XK15(3)  :  SPHI-SIN(PHI)  :  CPHI-COS(PHI) 

2430  NL(1)-XK15(2)^XK15(5)+QU*CTHETA*SPHI-A(1,3)*XK15t3) 
-2440  NL(2)-AIQRL*XK1 5(e)*XK1 5(4) 

2450  NL(3)-(XK1 5(6)*SPHI+XK1 5(4)*CPHI)‘TTHETA 
2460  NL(4)-AIPQN*XK15(2)*XK15(6) 

2470  NL(5)--XK15(2)‘XK15(1)*U/(U-A2ADOT)+32.17/(U-AZADOT) 
I  *(CTHETA*CPHI-COS(THETO)+SIN(THETO)‘XK(7)) 

_2480  NL(6)-AIRPM*XK15(4)*XK15(2)+AIPBM*NL(5) 

:~2490  NL(7)-XK15(6)*(CPHI-1)-XK15{4)*SPHI 
2500  RETURN 

2600  REM  S/R  FOR  STORING  DATA  TO  FILE 

2610  PRINT  #1,  USING  "###.#  i-###.#####  -I-###.##### 

2620  REM  CONTINUE  PRINT  USING  +###.#«#«#  +1t»UMUUU 
2630  REM  CONTINUE  +##».#####  •(•###.#####  •(-##«.#####"; 

2640  REM  CONTINUE  T.XK(1).XK(2).XK(3),XK(4),XK(5),XK(6).XK(7) 
2650  RETURN 

2700  REM  S/R  TO  CLOSE  DATA  FILE 
2710  CLOSE  #1 

2720  PRINT  #2,  USING  "+#.##  +1^##.#####  +###.##### 

2730  REM  CONTINUE  PRINT  USING  4###.#####";2eta,XKmax(1) 
2740  REM  CONTINUE  ,XKmax(3).XKmax(5) 

2750  CLOSE  #2 

2760  PRINT  "DATA  FILES  CLOSED" 

2770  RETURN 
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